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ABSTRACT 

We present time-dependent axisymmetric magnet ohydro dynamic simulations of 
the interaction of a relativistic magnetized wind produced by a proto-magnetar with a 
surrounding stellar envelope, in the first ~ 10 seconds after core collapse. We inject a 
super- magnetosonic wind with £^ = 10^^ ergs s~^ into a cavity created by an outgoing 
supernova shock. A strong toroidal magnetic field builds up in the bubble of plasma 
and magnetic field that is at first inertially confined by the progenitor star. This drives 
a jet out along the polar axis of the star, even though the star and the magnetar wind 
are each spherically symmetric. The jet has the properties needed to produce a long- 
duration gamma-ray burst (GRB). At ~ 5 s after core bounce, the jet has escaped the 
host star and the Lorentz factor of the material in the jet at large radii ~ 10^^ cm 
is similar to that in the magnetar wind near the source. Most of the spindown power 
of the central magnetar escapes via the relativistic jet. There are fluctuations in the 
Lorentz factor and energy flux in the j et on ~ 0.01 — 0.1 second timescale. These may 
contribute to variability in GRB emission (e.g., via internal shocks). 

Key words: Stars: neutron; stars: supernovae: general; gamma-rays: bursts; stars: 
winds, outflows; magnetic field; MHD 



1 INTRODUCTION 



Observations of long-duration gamma-ray bursts (GRBs) 
show that ma ny are assoc iated with core-collapse super- 
novae (SNe) (|Wooslev k BloomI l2006l). A newly formed 
rapidly rotating mag netar (e.g. JUsovlll992l : iThompsonll 1 994l : 
IWheeler et al. l200Qh or a black hole with an accretion disk 



fe.g.. iMacFad ven k Woosley| ll999l ) are possible central en- 
gines powering GRBs and their associated SNe. In both 
of these models, a key problem is to understand how the 
relativistic material necessary for generating the GRB es- 
capes from deep within the host star. A natural possibil- 
ity suggested by both afterglow observations ("jet breaks"; 
iRho ads 1999) and GRB energetics, is that a collimated out- 
flow punches through the stellar envelope, providing a chan- 
nel out of wh ich the relativistic material can escape (e.g., 
iMatznerlliooi ). 

In the collapsar model, collimated outflows from GRBs 
are accounted for by jets produced by an accretion flow onto 
a central black hole. In the magnetar model, the origin of 
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such collimated outflows is less clear because relativistic out- 
flows by themselves do no t e fficiently self-co llimate (e.g., 
iLvubarskv k EichleJ fiooih . iKonigl GranotI tm± sug- 
geste d, by analogy to pulsar wind nebulae (|Begelman k Lil 
I1992I ) (PWNe), that the interaction of the wind from the 
spinning-down magnetar with the surrounding star could fa- 
cilitate collimation. In a previous paper, we used axisymmet- 
ric thi n shell calculations to dem on strate that this can in fact 
occur (iBucciantini et af]|2007l : see lUzdenskv MacFadvenI 
I2OO6I . I2OO7I for related ideas based on force-free rather 
than inertially loaded outflows). As Konigl and Granot an- 
ticipated, the physical picture is analogous to that used 
to understand the morphologies of pulsar wind neb ulae 
(|Komissarov Lvubarsk^l2004l : [Pel Zanna et"alll2004l ): the 
magnetar wind shocks on the surrounding (exploding) stel- 
lar envelope, creating a bubble of relativistic plasma and 
magnetic field inside the host star (a "magnetar wind neb- 
ula"; MWN). If the toroidal magnetic field in the bubble 
is sufficiently strong, the bubble expands primarily in the 
polar direction due to the "tube of toothpaste" effect - the 
magnetic hoop stress in the equator creates a pressure gra- 
dient along the axis which is not balanced by any restoring 
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force, thus driving the flow preferentially in the polar direc- 
tion. The nebula itself is ultimately confined by the inertia 
of the SN envelope, to which little energy is transferred, in 
contrast to the pressure confinement in trad itional magnetic 
tower models (|Sherwin Lvnden-Bellll20Q7l l or the case of 
confi nement by a pressurized cocoo n inside the progenitor 
star (jUzdenskv MacFadvenll2006l ). 

Our previous results were based on a thin-shell model 
for the interaction between a magnetar wind and its host 
star. In this Letter, we build on this work by carrying 
out time- dependent axisymmetric magnetohydrodynamic 
(MHD) simulations of the interaction between a magnetar 
wind and a surrounding star. This paper is organized as a 
proof of principle to show that a mechanism developed in the 
PWNe context to explain the origin of observed jet/plume 
features, can also apply to proto-magnetars. For this reason 
we have selected a simplified set of parameters, deferring a 
detailed parameter study to future work. 

We assume that an outgoing SN shock has already cre- 
ated a central evacuated cavity and that the host star is 
spherically symmetric. A magnetar is assumed to inject en- 
ergy at a constant rate of ^ = 10^^ ergs s~^, a spindown 
power that is appropriate for proto-magnetars with P ^ 1 
ms and B - 10^^ G fe.g.. iMetzger et ahllioOTh . The inter- 
action between the magnetar wind and the host star de- 
pends on the strength of the toroidal magnetic field B in 
the MWN; this in turn depends on the magnetization in 
the wind at large radii, a = r^B^c/E. More specifically, 
cr{RTs) = (JTS, the magnetization at radii typical of the 
termination shock's radius Rts, is the essential control pa- 
rameter. At early times after core-bounce (< a few sec), 
young magnetars lose significant amounts of mass to a neu- 
trino driven wind, leading to ctlc <^ 1, where ctlc is the 
magnetization evaluated at the light cylinder. However, the 
mass loss rate decreases as the neutron star cools and the 
wind at later times 10 — 100 sec) is ultra-relativist ic 
with aLC ^ 10^-10^ ([Thompson et al.ll2004l : iMetzger etHI 
I2OO7I I. This late-time outflow is the most promising for pro- 
ducing a GRB. The relation between aLC - which can be 
calculated with some confidence - and a at larger radii in 
the wind is difficult to determine because of uncertainties in 
the conversion of magnetic energy into kin etic energy in rel- 
ativistic winds fsee lBucciantini et al.|[2007l for a detailed dis- 
cussion of this problem in the current context). In this Letter 
we focus on later times when aLC is large. We assume that 
non-ideal MHD processes analogous to those that operate in 
the PSR-PWN problem are also at work in proto-magnetar 
winds; these convert the high aLC wind at the light-cylinder 
into a moderate ars, high 7 wind at the free wind's termina- 
tion radius. Specifically, we consider relativistic winds with 
'jw = 10 or 25 and ars = 0.1; GRBs require somewhat 
higher 7 outflows but these are more difficult to simulate 
for numerical reasons. 

The remainder of this paper is organized as follows. In 
the next section (32|) we describe our numerical setup. We 
then present the results of our simulations (33|) and discuss 
their implications for models of GRBs (SI). 

2 NUMERICAL SETUP 

All of the simulations were performed using the shock- 
capturi ng central-scheme for rela t ivistic i deal MHD devel- 
oped Zanna Bucciantinil (|2002l ) fc lPef Zanna et al.l 



(|2003l ). using ideal gas EoS with adiabatic coefficint 4/3. We 
refer the reader to these papers for a detailed description of 
the equations and the numerical algorithms. 

The interaction of the magnetar wind with the sur- 
rounding SN progenitor is investigated by performing 2-D 
axisymmetric simulations on a spherical grid. The domain 
in is the first quadrant from ^ = 0to^ = 7r/2, with re- 
flecting boundary conditions for ve and B at the polar axis 
to enforce axisymmetry, and similar boundary conditions in 
the equatorial plane. The grid is uniform in the direction 
with 100 cells. Given that one needs to study a wide range of 
spatial scales, from the termination shock at ^ 10^ cm to the 
outer edge of the star at 2 10^° cm, we chose a grid that is 
logarithmic in radius and that extends from rmin = 10^ cm 
to Tmax = 5 10^*^ cm or Tmax = 10^^ cm, with 100 cells per 
decade in radius. Zeroth order extrapolation is assumed at 
the outer boundary. The code is second order in both space 
and time, with a monotonized central limiter, chosen in or- 
der to resolve the large density jump between the lighter 
relativistic plasma inside the MWN and the heavier stellar 
envelope (the density can increase by a factor of ^ 10^). 

For proto-magnetars with millisecond rotation peri- 
ods, the location of the Alfvenic surface is at ^ 10^ cm 
and the fast magn etosonic surface is at ~ 10^ — 10^ cm 
([Bucciantini et al] [2006: hereafter BOG). In the present sim- 
ulations, we do not attempt to resolve the dynamics of the 
wind in the sub-magnetosonic region close to the neutron 
star. Instead, at rmin we inject a super-magnetosonic wind 
with a fixed Lorentz factor oi — 10 or 7^o = 25 and a 
magnetization of a = 0.1, a value appropriate to distances 
on the order of the termination shock radius (the simulations 
were less stable for higher a); under these assumptions, a is 
conserved thoughoutb the upstream region. The wind is as- 
sumed to be cold with pc^ /p =100, and to contain a purely 
toroidal magnetic field. The energy flux E — 10^^ erg s~^ 
is kept constant during the entire simulation. For simplic- 
ity, we neglect latitudinal variations in the wind and assume 
that it is isotropic. Self-consistency requires that the termi- 
nation shock produced by the interaction of the magnetar 
wind with the surrounding star and MWN lie at > rmin- 
If not, the assumption of a super-magnetosonic wind would 
not be valid. We discuss this constraint more in the next 
section. 

W e use a SSM© stellar progenitor from IWooslev et aP 
(|2002l ). The wind is injected inside a cavity of radius ~ 10 
cm, which is roughly the size of the collapsing region in the 
first second after core bounce (our initial time). In order 
to simulate the effect of a SN shock propagating inside the 
progenitor, the region between 10^ and 210^ cm IS given 
an outward velocity corresponding to a total kinetic energy 
~ 2 10^^ erg, similar to that used in the one-dime nsional ex- 
plosion calculations of IWooslev Weaver! (|l995l ). The sur- 
face of the progenitor star is located at 2.5 10^° cm. We 
assume that the density outside the star falls off as r~^; for 
the range of radii we simulate, the results are independent 
of the outer density profile. 

3 RESULTS 

Figure [1] shows the density, pressure, flow velocity and flow 
streamlines for our simulations with — 10 at 3 differ- 
ent times, 4, 5, & 6 seconds after core bounce. At all times 
there is an axial high velocity jet subtending a few degrees 
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Figure 2. Inner velocity structure (in units of c) at 4 sec after 
core bounce. Distances are in 10^ cm. The termination shock is 
closer to the source along the pole because the pressure in the 
bubble is larger at the pole than at the equator. The flow pattern 
beyond the termination shock shows the formation of an axial 
jet /plume due to magnetic hoop stresses. Because the magnetic 
fleld goes to at the equator by symmetry, there is a small region 
the equator where material can escape in an equatorial channel. 
This has no effect on the evolution of the system. 



surrounded by a somewhat less relativistic cocoon subtend- 
ing 10°. The jet carries a significant fraction of the spin- 
down energy of the magnetar. Recall that the energy input 
in the simulation is spherically symmetric, as is the progeni- 
tor star. The aspherical evolution of the system is generated 
self-consistently by the dynamics of the MWN created inside 
the star. Note that the boundary between the MWN and the 
host star is particularly easy to identify via the jump in the 
pressure. 

To help explain the physical origin of the jet in Fig- 
ure [1] Figure [2] zooms in and shows the flow velocity and 
streamlines in the inner region (< 4 10^ cm) at 4 seconds. 
The high velocity v ^ c region at small radii in Figure [2] (in- 
dicated in red) is the freely expanding wind with — 10. 
The wind goes through a termination shock at '-^ 10^ cm due 
to the high bounding pressure of the MWN. Even though 
the energy flux in the freely expanding wind is isotropic, 
the termination shock itself is highly asymmetric, with the 
radius of the termination shock being significantly smaller 
along the pole. Physically, this is because the strong toroidal 
field in the MWN creates an anisotropic pressure distribu- 
tion with the press ure at the pole significantly larger than 
that at the equator (jBegelman Li|[l992l : lBucciantini et al.l 
I2OO7I V The higher polar pressure in turn causes the termi- 
nation shock to be at a smaller radius. For the simulations 
in Figures [1] & [21 with <j — 0.1 in the wind, we find that 
the equatorial termination shock always lies outside 10^ 
cm. Thus our injection of a super-magnetosonic wind is self- 
consistent. By contrast, in simulations with a larger value 
of cr in the wind, we found that the termination shock can 



move in to < rmin = 10^ cm. Because the termination shock 
would then lie within the Alfven surface, the spindown of the 
magnetar would be modified and our assumption of super- 
magnetosonic injection would be invalid. We defer a more 
detailed study of this interesting possibility to future work. 

In the post termination shock region. Figure [2] shows 
that the flow undergoes a large-scale circulation with much 
of the matter being diverted from the equatorial region to 
the pole where it flows out along the jet, just as in anal- 
ogo us calculations for the structure of pulsar wind nebu- 
lae feomissaro v Lvubarsk^l2004l : IPel Za nna et al.l l2004l ) . 
This is caused by the toroidal magnetic field that builds up 
in the MWN, which is nearly in e quipartition near the lo - 
cation where the back- flow starts (jPel Zanna et al.ll2004l l. 
Although collimation by hoop stresses is ineffective for the 
initially relativistic flow in the free wind, the non-relativist ic 
V ^ 0.5c flow in the post termination shock region is effec- 
tively collimated by the magnetic field. 

The MWN forms in about a sound crossing time ^ 
0.1 — 0.2 sec. It takes ^ 1 second for the magnetar to fill 
the MWN with enough energy for the magnetic stresses to 
become dynamically important, and about 1 more second 
for the jet/plume to emerge from the SN shock, still inside 
the star. The anisotropic pressure distribution in the MWN, 
with a much larger pressure at the pole than at the equa- 
tor, and the significant energy and momentum flux in the 
axial jet, act in concert to push out through the surround- 
ing star along the polar direction. At t = 4 sec, the MWN 
and jet are still fully contained within the star (Fig. [T]). By 
t = 5 sec, however, they have reached the radius of the star 
at ~ 2.5 10^*^ cm and at t = 6 seconds, the material in the 
jet at large radii has escaped the star and has accelerated 
back up to V ^ c. Note that the high velocity core of the jet 
is surrounded by a cocoon of less relativistic material that 
includes shocked stellar wind m aterial, as is generically ex- 
pected to be the case (jUzdenskv MacFadvenll2006l ). The 
opening angle of the high 7 core of the jet is a few degrees. 
Because the jet is only marginally resolved in our fiducial 
simulations, we ran a higher resolution simulation with twice 
the number of grid points in the theta direction near the axis 
of the jet; the results were nearly identical to those described 
here. 

Although the magnetic field is crucial for generating 
and colUmating the outflow seen in Figure [1] the field is not 
energetically dominant; the ratio of the magnetic energy to 
the thermal energy in the bubble is typically ~ 0.1 — 0.2, 
although it can reach ^ 1 in the region where the hoop stress 
is most effective (Fig. [2]). The enthalpy of the shocked gas in 
the post-termination-shock region primarily determines the 
asymptotic Lorentz factor in the jet. In turn, the enthalpy 
of the shocked gas is determined (through the termination 
shock's thermalization of the magnetar wind) by the wind 
Lorentz factor just upstream of the TS, 'jw Thus we find 
a roughly one-to-one relation between the Lorentz factor of 
the jet core at large radii and jw- for 7^ = 10 and rmax = 
5 X 10'^ cm 7 reaches ~ 5 in the jet, whereas for 'jyj — 25 
and Tmax = 10^^ cm, 7 15 at large radii. The upper panel 
in Figure [3] shows 7 on the axis of the jet as a function of 
radius for this simulation (at t = 9 sec), demonstrating the 
acceleration from the marginally relativistic plasma at small 
radii near the termination shock to the highly relativistic 
flow with 7 7^o at large radii. The acceleration is roughly 
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Figure 1. Evolution of a magnetized bubble inflated by a magnetar wind with E = lO^-'^ergs s"-*^, = 10, and a = 0.1, inside a SSM© 
progenitor star. From left to right: density (g cm~^), pressure (g cm~^ c^) and velocity (in units of c). From top to bottom: snapshots at 
4, 5, & 6 seconds after core bounce. Distances are in 10^ cm; the radius of the progenitor star is 2.5 lO-*^*-* cm. By t = 5 s (middle panel) 
the jet has escaped the progenitor star. 



linear in radius, consistent with that expected in ; 
with a constant opening angle. 



"fireball" that a significant fraction of the energy is carried by the high 
velocity core of the jet. 



The lower panel in Figure [3] shows the ratio of the power 
in the high velocity core of the jet (solid) and in the wider 
cocoon (dashed) relative to the power supplied by the central 
source E; as in the upper panel, these calculations are for 
7to = 25 at t = 9 sec when the outflow has escaped the star. 
Figure [3] shows that at late times nearly all of the energy 
supplied by the central magnetar escapes to large radii and 



The upper and lower panels of Figure [3] demonstrate 
radial fluctuations in the Lorentz factor and E, respectively, 
despite the fact that the central source is assumed to have 
an E and that are constant in time. These fluctuations 
arise largely from shear instability near the base of the jet, 
and changes in the circulation of matter in the MWN (the 
region seen in Fig. [2]) - this accounts for the characteristic 
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R(1 0^ cm) 

Figure 3. Upper panel: Lorentz factor on the axis of the jet, 9 sec 
after core bounce for = 25 and rmax = lO-*^-*^ cm. 7 increases 
roughly Unearly in radius. Lower panel: Ratio of the radial energy 
flux in the high velocity core of the jet to the total source power 
(solid line). Also shown is the ratio of the total radial energy flux 
in the cocoon with respect to the source power (dashed line). 



length scale of ^ 10^ cm in Figure [3] The corresponding 
variability timescale is ^ 0.01 — 0.1 sec. 

4 DISCUSSION & IMPLICATIONS 

Our simulations demonstrate that a wind with properties 
typical of proto- magnetar winds inside a massive star can 
be coUimated into an axial jet, even in the conservative case 
in which the host star and wind are spherically symmetric. 
Physically, this occurs because of the strong toroidal mag- 
netic field in the bubble of relativistic plasma and magnetic 
field that is at first inertially confined by the progenitor star 
and SN shock. The magnetic field creates an asymmetric 
pressure distribution in the bubble, with the pressure much 
larger at the pole than the equator. In addition, the toroidal 
field collimates the moderately relativistic flow behind the 
termination shock into an axial jet. 

Although we have not been able to carry out simulations 
with Lorentz factors of ~ 100 — 10^ (largely for numerical 
reasons), as is required to explain GRBs, extrapolating our 
results suggests that an isotropic or equatorial flow with 
7 ~ 100 — 10^ would produce a jet at large radii with a 
comparable Lorentz factor. This collimated outflow would 
thus have physical properties similar to those required to 
produce GRBs. As noted in §1, calculations of the mass- loss 
rate from newly formed magnetars flnd that magnetar winds 
naturally hav e cflc > 100—10^ roughly 10-100 seconds af ter 
core bounce ([Thompson et al .^ 2004; Me tzger et"alll2007l l. 

Figure [3] shows that the jet escaping the star at large 
distances has fluctuations in 7 and ^ on ^ 0.01 — 0.1 sec- 
ond timescales. It is natural to speculate that this variabil- 
ity in the source region could manifest itself as variabil- 
ity in the gamma-ray emission in GRBs (via, e.g., inter- 
nal shocks). This is in addition to any intrinsic variability 
in the spindown of the central magnetar. It is also worth 
noting that the magnetic fleld in the MWN will be ad- 
vected out to la rge radii where the GRB emission occurs 
(jLyutikov Blandford 2003); the resulting fleld is likely 
strong enough to account for the observed gamma-ray emis- 
sion via synchrotron radiation, without the need for shock 
generated magnetic flelds. 

This Letter has focused on one part of the param- 



eter space of proto- magnetar and host star interactions, 
namely a magnetar wind with cr 0.1 spinning down into 
a cavity evacuated by a successful core- collapse SN. We 
flnd that in simulations with a more strongly magnetized 
wind (cr > 0.1), the termination shock moves inside the 
Alfvenic and fast magnetosonic surfaces (SSJ. The spindown 
of the central magnetar will thus be modifled. It is clearly 
of considerable interest to understand the coupled dynam- 
ics of the magnetar wind and the surrounding nebula in 
this limit. In addition, in our current simulations, nearly 
all of the spindown power of the central magnetar is chan- 
neled through the polar jet (Fig. [3]). Thus the late-time spin- 
down of the magnetar can Ukely generate a GRB, but it 
will not energize the surrounding SN. It is thus not clear 
if the current model produ ces both a hyper- energet ic SN 
and a GRB, as is observed ( Wooslev BloomI 120061 ) . One 
possibility is that the initial explosion itself is highly ener- 
getic because rotational ener gy contributes to the explosion 
fe. g., iThompson et al. 20051 ). Alternatively, we speculated 



in iBucciantini et al.l ([20071) that because the initial magne- 
tar wind is likely to be signiflcantly mass loaded and non- 
relativistic, it may contribute to energizing the SN shock 
rather than to producing a GRB. This possibility will be 
investigated in future work with simulations that properly 
take into account the evolution of E and a with time in mag- 
netar winds. Finally, it is clearly important to understand 
the stability of the strong toroidal fleld that is crucial to the 
dynamics of the MWN and the axial jet. This can be studied 
both analytically and with three dimensional simulations. 

Near the co mplet ion of this work, 

Komissarov Bar kovl (|2007[ ) submitted a paper on 
magnetar spindown inside a star whose results are com- 
plementary to ours. They calculated the dynamics of an 
outflow from a magnetar interacting with an infalling 
stellar progenitor for 200 ms after core-collapse, flnding a 
collimated non-relativistic jet-like outflow. By contrast, we 
inject a super-magnetosonic wind with properties derived 
from separate spindown calculations, and are thus able to 
simulate the evolution for a much longer period of time ~ 10 
sec and to much higher Lorentz factors. Although a number 
of the details differ, the broad conclusions from these two 
works are similar: magnetar formation can produce a jet 
with properties similar to those required to produce GRBs. 
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